%%Test Program for Strang Voltage Divider


% http://cns.bu.edu/~eric/readings/strang.pdf

%% Circuit Diagram
% http://cns.bu.edu/~eric/readings/VD.jpg
% 


%%  Variables
% The following matrixes and vectors are needed
%
% Sparse matrixes are used throughout the code
%
% * A Incidence matrix  node by edge
% * A_0 Reduced incidence matrix nodes by edge-1
% * b Voltage source vector b edge by 1
% * f Current source vector node by 1
% * G Conductance matrix G edge by edge
% * C Capacitance matrix edge by edge
% * y Current vector  node by 1
% * x Potentials node by 1
% * e Potential differences edge by 1

%% Kirchoff KCL and KVL
% KCL and KVL may be succinctly stated in matrix vector form as  follows.
%
%% Ohms law (constitutive) appears below in conductance form: 
% $$y_r=Ge$$
%

%% A capacitive constituitive law of the form
%
% $$y_c=\mbox{-}C\dot{e}=\mbox{-}C\frac{d}{dt}(b\mbox{-}A_0\dot{x})~=~C\dot{e}$$
%


%% KCL with substitutions of KVL, Ohmic Current and Displacement Current
%
% $$A_0^tG(b~ \mbox{-}~A_0x)+A_0^t CA_0\dot{x} ~=~ f$$
%
%% SUMMARY OF KCL/KVL/OHM
%% Voltage Drop on Edges
% 
% $$e=b~ \mbox{-}~ A_0x$$
%

%% Ohms Law for resistive currents
% $$y_r=Ge$$

%% KCL balances to current sources f
%
% $$A_0^ty=f$$

%% KCL/KVL with only resistive edge variables according to Strang
%%
% $$A_0^tG(b~\mbox{-}~A_0x)=f$$
%

% ELS 11/1/2006


%%Incidence matrix
% A= -1 1;1 -1
% A_0= [-1 1]^t
i=[];j=[];s=[]; %initialize to empty
%
i=[i 1]; j=[j 1]; s=[s -1];
i=[i 1]; j=[j 2]; s=[s 1];
i=[i 2]; j=[j 1]; s=[s 1];
i=[i 2]; j=[j 2]; s=[s -1];
A=sparse(i,j,s);
%
%% Reduced incidence matrix: chop last column
A_0=A(:,1); % chop last column to form reduced Incidence matrix
%
%% Conductivity matrix
R1=10; %10 Ohms
R2=5; %5 Ohms
G=[1/R1 0;0 1/R2];
G=sparse(G);
%% Voltage sources
%Battery V_0 pos up
%two Resistors R_1=10 and R_2=5
%
V_0=10; % 10 volt source
b=[-V_0;0]; % b is a column vector
b=sparse(b);
%% Solve

%%
%
% $$A_0^t G (b~ \mbox{-}~A_0 x) ~=~ 0$$
% 

%%
%
% $$x~=~({A_0}^t G A_0)^{\mbox{-1}}A_0^t G b$$
%
x=A_0'*G*b/(A_0'*G*A_0)
